In this study, based on DeepONet [1] a Physics-Informed Neural Networks (PINNs), which has been validated in solving both ordinary differential equations (ODEs) and partial differential equations (PDEs), we aim to train a neural network that represents the mapping from the inputs, including horizontal coordinate x, vertical coordinate z and time t, to the corresponding outputs, which include potential ϕ and elevation of the free surface η. In this architecture, the Laplace equation as a governing equation, along with the boundary conditions, especially, the dynamic and kinematic equations, forms part of the loss function to ensure that the results comply with the laws of physics. The training dataset is obtained through the numerical results of the Numerical Wave Tank model (NWT) [2]. During the training of Neural Networks, especially PINNs, it is often challenging to ensure both accuracy and convergence. One of the difficulties encountered is the spectral bias. Due to the gradient descent optimization, neural networks tend to prioritize learning low-frequency information (as low-frequency components have a greater impact on gradients during backpropagation), while neglecting high-frequency information. To resolve this difficulty, Wang et al. 2021[3] has proposed the Fourier Feature method, which construct a coordinate embedding by using an original formulation. In wave propagation, the high-frequency information, e.g., the wave breaking, play an important role, therefore, Fourier Features method is embedded as the preprocessing step of the input data into the DeepONet, and the results are compared with those obtained from fully connected neural networks.